Geospatial/SAE workshop
Small area estimation
Josh Merfeld
November 2025
SAE in Malawi
Let’s take a look at Malawi
Consider the 2019/2020 Integrated Household Survey (IHS5)
- Was used for the Malawi Poverty Report 2020
- Can say things about poverty at the district level
- If you want to split by urban/rural, only at the region level
- So what if we want to estimate poverty below the district, at the TA (admin3)?
Three “flavors” of SAE
Estimate a model at the household level, aggregate to admin3
- This is a “unit-level” model
Estimate a motel at the admin4 level, aggregate to admin3
- This is also a “unit-level” model, but some people strongly disagree and prefer to call this a “sub-area” model
Estimate a motel at the admin3 level
- This is an area-level model
Household-level model
- Traditionally at the World Bank, SAE was done at the household level
- Elbers, Lanjouw, and Lanjouw (ELL) method
- ELL is a completely “synthetic” estimator
- (Can also use EBP for this, more in a bit)
- Requirements:
- Household survey data with predictors (e.g. asset ownership) and outcome (e.g. expenditures/poverty)
- Household census data with the same predictors
- How many of your countries have this?
“Sub-area” model
- Option two is to estimate a model at the admin4 level
- Aggregate to admin3 level
- ELL could in theory work here but it’s not usually done
- EBP most common in our work
- Requirements:
- Survey outcome at the admin4 level; this can be from household survey data or some other type of survey
- Predictors at the admin 4 level; we can use geospatial data here!
EBP (empirical best predictor)
- “Subarea model” (variant of Battese, Harter, Fuller 1988)
\[g(\hat{p}_{ic}^{dir})=X_{ic}'\beta+b_i+v_i+\varepsilon_{ic}\]
- \(\hat{p}_{ic}^{dir}\) is poverty at the cluster (subarea) level, from a survey
- \(g(\cdot)\) is some transformation of poverty (we prefer arcsin transformation)
- \(X_{ic}'\) is a vector of predictors
- \(v_i\sim N(0,\sigma_v^2)\) and \(\varepsilon_{ic}\sim N(0,\sigma_\varepsilon^2)\)
- \(\beta\), \(\sigma_v^2\), and \(\sigma_\varepsilon^2\) are estimated parameters from a mixed model
EBP (empirical best predictor)
- Key thing: EBP combines the survey estimates with predictions from the model
- Intuition: areas that have larger sample size will put more weight on the sample than the prediction
- No sample in an area? Synthetic prediction only!
- Uncertainty estimated through a parametric bootstrap
- Won’t go into detail here, but…
- Assumes normality of the random effect and the error term
- Research suggests it is robust to reasonable deviations from normality
“Area-level” model
- Most common: Fay-Herriot
\[y_i = X_i'\beta+v_i+e_i\]
- \(X_i\) vector of predictors
- \(v_i\) model-based random effect
- \(e_i\) design-based random errors
- Problem:
- We usually need to give the model the variance of the outcome
- We use survey for this, but some areas don’t have survey observations!
- Some areas don’t have variance in the outcome, even with observations!
Feature selection
- Let’s start with some example data I have
- This comes from Malawi
- Northern Malawi only (due to the size of the data)
- And we’ll use it all day tomorrow!
- This comes from Malawi
A short explanation
- The survey data is collapsed to the admin3 level (TAs)
- This is the
area, in SAE terminology - I have poverty rates for
areas(TAs) andsubareas(EAs) - I have some variables that predict poverty at the
subarealevel
- This is the
- So it’s a perfect setup for SAE!
- We want to estimate poverty at the TA
- We don’t have any observations in some TAs and we have too few in others
- We could estimate a subarea model
Observations?
Predictive features
- I also have a bunch of predictive features!
- The data come from something called MOSAIKS, that we’ll discuss briefly tomorrow
- In short, they are variables derived from satellite imagery
- Take a look at this
# A tibble: 2,911 × 501
EA_CODE mosaik1 mosaik2 mosaik3 mosaik4 mosaik5 mosaik6 mosaik7 mosaik8
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 10101001 0.00143 0.00242 0.632 0.0334 0.0684 0.223 0.00641 0.0753
2 10101002 0.000659 0.000250 0.911 0.0468 0.0978 0.357 0.00359 0.127
3 10101003 0.000657 0.000403 0.811 0.0373 0.0794 0.326 0.00285 0.100
4 10101004 0.00102 0.000769 0.975 0.0578 0.111 0.369 0.00584 0.140
5 10101005 0.000472 0.000351 0.815 0.0344 0.0668 0.381 0.00287 0.0947
6 10101006 0.00107 0.000835 0.861 0.0496 0.122 0.315 0.00536 0.137
7 10101007 0.00132 0.000842 1.13 0.0549 0.0999 0.594 0.00649 0.154
8 10101008 0.00202 0.00182 1.05 0.0796 0.166 0.415 0.00953 0.179
9 10101009 0.000445 0.000417 0.834 0.0332 0.0663 0.375 0.00278 0.0950
10 10101010 0.000720 0.000438 0.794 0.0367 0.0849 0.328 0.00377 0.109
# ℹ 2,901 more rows
# ℹ 492 more variables: mosaik9 <dbl>, mosaik10 <dbl>, mosaik11 <dbl>,
# mosaik12 <dbl>, mosaik13 <dbl>, mosaik14 <dbl>, mosaik15 <dbl>,
# mosaik16 <dbl>, mosaik17 <dbl>, mosaik18 <dbl>, mosaik19 <dbl>,
# mosaik20 <dbl>, mosaik21 <dbl>, mosaik22 <dbl>, mosaik23 <dbl>,
# mosaik24 <dbl>, mosaik25 <dbl>, mosaik26 <dbl>, mosaik27 <dbl>,
# mosaik28 <dbl>, mosaik29 <dbl>, mosaik30 <dbl>, mosaik31 <dbl>, …We have a problem
- What’s the problem?
- It’s actually impossible to estimate a model with more predictors than observations!
Another problem: overfitting
There’s another problem, too
If we have too many predictors, we can “overfit” the model
- This means the model is too complex
- It fits the data we have too well
- This means it doesn’t generalize well to new data
So we need to select the best predictors
- What does “best” mean here?
Generalizing out-of-sample
We want to know what best predicts OUT of sample
So we are going to set up our data to allow this:
- We will split the data into X parts
- A common number for X is 10, but let’s do 5
Cross validation
Cross validation
Cross validation - random folds
Code
# A tibble: 6 × 5
EA_CODE poor total_weights total_obs fold
<dbl> <dbl> <dbl> <dbl> <int>
1 10101006 0.230 5690. 16 3
2 10101011 0.444 7614. 16 5
3 10101027 0.0947 9441. 16 4
4 10101033 0.376 7486. 16 3
5 10101039 0.600 9147. 16 4
6 10101054 0.497 5351. 16 5Cross validation
Cross validation
But what “models” are we going to fit?
- What are the models we are going to fit?
- We want a way to select the best predictors
- This will reduce the number of predictors and prevent overfitting (we hope)
- We are going to use a method called LASSO (or lasso)
- It’s an acronym: Least Absolute Shrinkage and Selection Operator
- No details, but it’s a way to select the best predictors
- It “penalizes” the coefficients of the predictors
Rpackageglmnetdoes this for us
The setup - with a transformed outcome
Code
library(glmnet)
set.seed(398465) # this is a random process, so we want to set the seed!
# we need to set up the data (combining the predictors and the outcome)
data <- surveycollapsed |>
left_join(predictors, by = "EA_CODE")
# cv.glmnet will set up everything for us
lasso <- cv.glmnet(
y = asin(sqrt(data$poor)), # the outcome
x = data |> dplyr::select(starts_with("mosaik")) |> as.matrix(), # the predictors (as.matrix() is required)
weights = data$total_weights, # the weights (sample weights)
nfolds = 5) # number of folds (10 is the default)
lasso
Call: cv.glmnet(x = as.matrix(dplyr::select(data, starts_with("mosaik"))), y = asin(sqrt(data$poor)), weights = data$total_weights, nfolds = 5)
Measure: Mean-Squared Error
Lambda Index Measure SE Nonzero
min 0.02030 26 0.04227 0.006409 6
1se 0.06493 1 0.04418 0.005811 0What have we done?
Call: cv.glmnet(x = as.matrix(dplyr::select(data, starts_with("mosaik"))), y = asin(sqrt(data$poor)), weights = data$total_weights, nfolds = 5)
Measure: Mean-Squared Error
Lambda Index Measure SE Nonzero
min 0.02030 26 0.04227 0.006409 6
1se 0.06493 1 0.04418 0.005811 0- What are the different “models”?
- Different values of lambda
- In this case, the “best” lambda is 0.02030
- Note that some people prefer to use the
1sevalue (it is more conservative). No details today.
Different values of lambda: different predictors!
Call: cv.glmnet(x = as.matrix(dplyr::select(data, starts_with("mosaik"))), y = asin(sqrt(data$poor)), weights = data$total_weights, nfolds = 5)
Measure: Mean-Squared Error
Lambda Index Measure SE Nonzero
min 0.02030 26 0.04227 0.006409 6
1se 0.06493 1 0.04418 0.005811 0- At the “optimal” lambda, we have 6 predictors (non-zero coefficients)
Choosing based on mean-squared error (MSE)
Non-zero coefficients
Non-zero coefficients
501 x 1 sparse Matrix of class "dgCMatrix"
s1
(Intercept) 0.568658061
mosaik1 .
mosaik2 .
mosaik3 .
mosaik4 .
mosaik5 .
mosaik6 .
mosaik7 .
mosaik8 .
mosaik9 .
mosaik10 .
mosaik11 .
mosaik12 .
mosaik13 .
mosaik14 .
mosaik15 .
mosaik16 .
mosaik17 .
mosaik18 .
mosaik19 .
mosaik20 .
mosaik21 .
mosaik22 .
mosaik23 .
mosaik24 .
mosaik25 .
mosaik26 .
mosaik27 .
mosaik28 .
mosaik29 .
mosaik30 .
mosaik31 .
mosaik32 .
mosaik33 .
mosaik34 .
mosaik35 .
mosaik36 .
mosaik37 .
mosaik38 .
mosaik39 -0.003056014
mosaik40 .
mosaik41 .
mosaik42 .
mosaik43 .
mosaik44 .
mosaik45 .
mosaik46 .
mosaik47 .
mosaik48 .
mosaik49 .
mosaik50 .
mosaik51 .
mosaik52 .
mosaik53 .
mosaik54 .
mosaik55 .
mosaik56 .
mosaik57 .
mosaik58 .
mosaik59 .
mosaik60 .
mosaik61 .
mosaik62 .
mosaik63 .
mosaik64 .
mosaik65 .
mosaik66 .
mosaik67 .
mosaik68 .
mosaik69 .
mosaik70 .
mosaik71 .
mosaik72 .
mosaik73 .
mosaik74 .
mosaik75 .
mosaik76 .
mosaik77 .
mosaik78 .
mosaik79 .
mosaik80 .
mosaik81 .
mosaik82 .
mosaik83 .
mosaik84 .
mosaik85 .
mosaik86 .
mosaik87 .
mosaik88 .
mosaik89 .
mosaik90 .
mosaik91 .
mosaik92 .
mosaik93 .
mosaik94 .
mosaik95 .
mosaik96 .
mosaik97 .
mosaik98 .
mosaik99 .
mosaik100 .
mosaik101 .
mosaik102 .
mosaik103 .
mosaik104 .
mosaik105 .
mosaik106 .
mosaik107 .
mosaik108 .
mosaik109 .
mosaik110 .
mosaik111 .
mosaik112 .
mosaik113 .
mosaik114 .
mosaik115 .
mosaik116 .
mosaik117 .
mosaik118 .
mosaik119 .
mosaik120 .
mosaik121 .
mosaik122 .
mosaik123 .
mosaik124 .
mosaik125 .
mosaik126 .
mosaik127 .
mosaik128 .
mosaik129 .
mosaik130 .
mosaik131 .
mosaik132 .
mosaik133 .
mosaik134 .
mosaik135 .
mosaik136 .
mosaik137 .
mosaik138 .
mosaik139 .
mosaik140 .
mosaik141 .
mosaik142 .
mosaik143 .
mosaik144 .
mosaik145 .
mosaik146 .
mosaik147 .
mosaik148 .
mosaik149 .
mosaik150 .
mosaik151 .
mosaik152 .
mosaik153 .
mosaik154 .
mosaik155 .
mosaik156 .
mosaik157 .
mosaik158 .
mosaik159 .
mosaik160 .
mosaik161 .
mosaik162 .
mosaik163 .
mosaik164 .
mosaik165 .
mosaik166 .
mosaik167 .
mosaik168 .
mosaik169 .
mosaik170 .
mosaik171 .
mosaik172 .
mosaik173 .
mosaik174 .
mosaik175 .
mosaik176 .
mosaik177 .
mosaik178 .
mosaik179 .
mosaik180 .
mosaik181 .
mosaik182 .
mosaik183 .
mosaik184 .
mosaik185 .
mosaik186 .
mosaik187 .
mosaik188 .
mosaik189 .
mosaik190 .
mosaik191 .
mosaik192 .
mosaik193 .
mosaik194 .
mosaik195 .
mosaik196 .
mosaik197 .
mosaik198 .
mosaik199 .
mosaik200 .
mosaik201 .
mosaik202 .
mosaik203 .
mosaik204 .
mosaik205 .
mosaik206 .
mosaik207 .
mosaik208 .
mosaik209 .
mosaik210 .
mosaik211 .
mosaik212 .
mosaik213 .
mosaik214 .
mosaik215 .
mosaik216 .
mosaik217 .
mosaik218 .
mosaik219 .
mosaik220 .
mosaik221 .
mosaik222 .
mosaik223 .
mosaik224 .
mosaik225 .
mosaik226 .
mosaik227 .
mosaik228 .
mosaik229 .
mosaik230 .
mosaik231 .
mosaik232 .
mosaik233 .
mosaik234 0.015558794
mosaik235 .
mosaik236 .
mosaik237 .
mosaik238 .
mosaik239 .
mosaik240 .
mosaik241 .
mosaik242 .
mosaik243 .
mosaik244 .
mosaik245 .
mosaik246 .
mosaik247 .
mosaik248 .
mosaik249 .
mosaik250 .
mosaik251 .
mosaik252 .
mosaik253 .
mosaik254 .
mosaik255 .
mosaik256 .
mosaik257 .
mosaik258 .
mosaik259 .
mosaik260 .
mosaik261 .
mosaik262 .
mosaik263 .
mosaik264 .
mosaik265 .
mosaik266 .
mosaik267 .
mosaik268 .
mosaik269 .
mosaik270 .
mosaik271 .
mosaik272 .
mosaik273 .
mosaik274 .
mosaik275 .
mosaik276 .
mosaik277 0.044130256
mosaik278 .
mosaik279 .
mosaik280 -8.148026851
mosaik281 .
mosaik282 .
mosaik283 .
mosaik284 .
mosaik285 .
mosaik286 .
mosaik287 .
mosaik288 .
mosaik289 .
mosaik290 .
mosaik291 .
mosaik292 .
mosaik293 .
mosaik294 .
mosaik295 .
mosaik296 .
mosaik297 .
mosaik298 .
mosaik299 .
mosaik300 .
mosaik301 .
mosaik302 .
mosaik303 .
mosaik304 .
mosaik305 .
mosaik306 .
mosaik307 .
mosaik308 .
mosaik309 .
mosaik310 .
mosaik311 .
mosaik312 .
mosaik313 .
mosaik314 .
mosaik315 .
mosaik316 .
mosaik317 .
mosaik318 .
mosaik319 .
mosaik320 .
mosaik321 .
mosaik322 .
mosaik323 .
mosaik324 .
mosaik325 .
mosaik326 .
mosaik327 .
mosaik328 .
mosaik329 .
mosaik330 .
mosaik331 .
mosaik332 .
mosaik333 .
mosaik334 .
mosaik335 .
mosaik336 .
mosaik337 .
mosaik338 .
mosaik339 .
mosaik340 .
mosaik341 .
mosaik342 .
mosaik343 .
mosaik344 .
mosaik345 .
mosaik346 .
mosaik347 .
mosaik348 .
mosaik349 .
mosaik350 .
mosaik351 .
mosaik352 .
mosaik353 .
mosaik354 .
mosaik355 .
mosaik356 .
mosaik357 .
mosaik358 .
mosaik359 .
mosaik360 .
mosaik361 .
mosaik362 .
mosaik363 .
mosaik364 .
mosaik365 .
mosaik366 .
mosaik367 .
mosaik368 .
mosaik369 .
mosaik370 .
mosaik371 .
mosaik372 .
mosaik373 .
mosaik374 .
mosaik375 .
mosaik376 .
mosaik377 .
mosaik378 .
mosaik379 .
mosaik380 .
mosaik381 .
mosaik382 .
mosaik383 .
mosaik384 .
mosaik385 .
mosaik386 .
mosaik387 .
mosaik388 .
mosaik389 .
mosaik390 .
mosaik391 .
mosaik392 .
mosaik393 .
mosaik394 .
mosaik395 .
mosaik396 -35.261808063
mosaik397 .
mosaik398 .
mosaik399 .
mosaik400 .
mosaik401 .
mosaik402 .
mosaik403 .
mosaik404 .
mosaik405 .
mosaik406 .
mosaik407 .
mosaik408 .
mosaik409 .
mosaik410 .
mosaik411 .
mosaik412 .
mosaik413 .
mosaik414 .
mosaik415 .
mosaik416 .
mosaik417 .
mosaik418 .
mosaik419 .
mosaik420 .
mosaik421 .
mosaik422 .
mosaik423 .
mosaik424 .
mosaik425 .
mosaik426 .
mosaik427 .
mosaik428 .
mosaik429 .
mosaik430 .
mosaik431 .
mosaik432 .
mosaik433 .
mosaik434 .
mosaik435 .
mosaik436 .
mosaik437 .
mosaik438 .
mosaik439 .
mosaik440 .
mosaik441 .
mosaik442 .
mosaik443 .
mosaik444 .
mosaik445 .
mosaik446 .
mosaik447 .
mosaik448 .
mosaik449 .
mosaik450 .
mosaik451 .
mosaik452 .
mosaik453 .
mosaik454 .
mosaik455 .
mosaik456 .
mosaik457 .
mosaik458 .
mosaik459 55.635488670
mosaik460 .
mosaik461 .
mosaik462 .
mosaik463 .
mosaik464 .
mosaik465 .
mosaik466 .
mosaik467 .
mosaik468 .
mosaik469 .
mosaik470 .
mosaik471 .
mosaik472 .
mosaik473 .
mosaik474 .
mosaik475 .
mosaik476 .
mosaik477 .
mosaik478 .
mosaik479 .
mosaik480 .
mosaik481 .
mosaik482 .
mosaik483 .
mosaik484 .
mosaik485 .
mosaik486 .
mosaik487 .
mosaik488 .
mosaik489 .
mosaik490 .
mosaik491 .
mosaik492 .
mosaik493 .
mosaik494 .
mosaik495 .
mosaik496 .
mosaik497 .
mosaik498 .
mosaik499 .
mosaik500 . What we want: the non-zero variable names!
- Getting the names of the variables is more complicated than it should be
Code
s1
(Intercept) 0.568658061
mosaik1 0.000000000
mosaik2 0.000000000
mosaik3 0.000000000
mosaik4 0.000000000
mosaik5 0.000000000
mosaik6 0.000000000
mosaik7 0.000000000
mosaik8 0.000000000
mosaik9 0.000000000
mosaik10 0.000000000
mosaik11 0.000000000
mosaik12 0.000000000
mosaik13 0.000000000
mosaik14 0.000000000
mosaik15 0.000000000
mosaik16 0.000000000
mosaik17 0.000000000
mosaik18 0.000000000
mosaik19 0.000000000
mosaik20 0.000000000
mosaik21 0.000000000
mosaik22 0.000000000
mosaik23 0.000000000
mosaik24 0.000000000
mosaik25 0.000000000
mosaik26 0.000000000
mosaik27 0.000000000
mosaik28 0.000000000
mosaik29 0.000000000
mosaik30 0.000000000
mosaik31 0.000000000
mosaik32 0.000000000
mosaik33 0.000000000
mosaik34 0.000000000
mosaik35 0.000000000
mosaik36 0.000000000
mosaik37 0.000000000
mosaik38 0.000000000
mosaik39 -0.003056014
mosaik40 0.000000000
mosaik41 0.000000000
mosaik42 0.000000000
mosaik43 0.000000000
mosaik44 0.000000000
mosaik45 0.000000000
mosaik46 0.000000000
mosaik47 0.000000000
mosaik48 0.000000000
mosaik49 0.000000000
mosaik50 0.000000000
mosaik51 0.000000000
mosaik52 0.000000000
mosaik53 0.000000000
mosaik54 0.000000000
mosaik55 0.000000000
mosaik56 0.000000000
mosaik57 0.000000000
mosaik58 0.000000000
mosaik59 0.000000000
mosaik60 0.000000000
mosaik61 0.000000000
mosaik62 0.000000000
mosaik63 0.000000000
mosaik64 0.000000000
mosaik65 0.000000000
mosaik66 0.000000000
mosaik67 0.000000000
mosaik68 0.000000000
mosaik69 0.000000000
mosaik70 0.000000000
mosaik71 0.000000000
mosaik72 0.000000000
mosaik73 0.000000000
mosaik74 0.000000000
mosaik75 0.000000000
mosaik76 0.000000000
mosaik77 0.000000000
mosaik78 0.000000000
mosaik79 0.000000000
mosaik80 0.000000000
mosaik81 0.000000000
mosaik82 0.000000000
mosaik83 0.000000000
mosaik84 0.000000000
mosaik85 0.000000000
mosaik86 0.000000000
mosaik87 0.000000000
mosaik88 0.000000000
mosaik89 0.000000000
mosaik90 0.000000000
mosaik91 0.000000000
mosaik92 0.000000000
mosaik93 0.000000000
mosaik94 0.000000000
mosaik95 0.000000000
mosaik96 0.000000000
mosaik97 0.000000000
mosaik98 0.000000000
mosaik99 0.000000000
mosaik100 0.000000000
mosaik101 0.000000000
mosaik102 0.000000000
mosaik103 0.000000000
mosaik104 0.000000000
mosaik105 0.000000000
mosaik106 0.000000000
mosaik107 0.000000000
mosaik108 0.000000000
mosaik109 0.000000000
mosaik110 0.000000000
mosaik111 0.000000000
mosaik112 0.000000000
mosaik113 0.000000000
mosaik114 0.000000000
mosaik115 0.000000000
mosaik116 0.000000000
mosaik117 0.000000000
mosaik118 0.000000000
mosaik119 0.000000000
mosaik120 0.000000000
mosaik121 0.000000000
mosaik122 0.000000000
mosaik123 0.000000000
mosaik124 0.000000000
mosaik125 0.000000000
mosaik126 0.000000000
mosaik127 0.000000000
mosaik128 0.000000000
mosaik129 0.000000000
mosaik130 0.000000000
mosaik131 0.000000000
mosaik132 0.000000000
mosaik133 0.000000000
mosaik134 0.000000000
mosaik135 0.000000000
mosaik136 0.000000000
mosaik137 0.000000000
mosaik138 0.000000000
mosaik139 0.000000000
mosaik140 0.000000000
mosaik141 0.000000000
mosaik142 0.000000000
mosaik143 0.000000000
mosaik144 0.000000000
mosaik145 0.000000000
mosaik146 0.000000000
mosaik147 0.000000000
mosaik148 0.000000000
mosaik149 0.000000000
mosaik150 0.000000000
mosaik151 0.000000000
mosaik152 0.000000000
mosaik153 0.000000000
mosaik154 0.000000000
mosaik155 0.000000000
mosaik156 0.000000000
mosaik157 0.000000000
mosaik158 0.000000000
mosaik159 0.000000000
mosaik160 0.000000000
mosaik161 0.000000000
mosaik162 0.000000000
mosaik163 0.000000000
mosaik164 0.000000000
mosaik165 0.000000000
mosaik166 0.000000000
mosaik167 0.000000000
mosaik168 0.000000000
mosaik169 0.000000000
mosaik170 0.000000000
mosaik171 0.000000000
mosaik172 0.000000000
mosaik173 0.000000000
mosaik174 0.000000000
mosaik175 0.000000000
mosaik176 0.000000000
mosaik177 0.000000000
mosaik178 0.000000000
mosaik179 0.000000000
mosaik180 0.000000000
mosaik181 0.000000000
mosaik182 0.000000000
mosaik183 0.000000000
mosaik184 0.000000000
mosaik185 0.000000000
mosaik186 0.000000000
mosaik187 0.000000000
mosaik188 0.000000000
mosaik189 0.000000000
mosaik190 0.000000000
mosaik191 0.000000000
mosaik192 0.000000000
mosaik193 0.000000000
mosaik194 0.000000000
mosaik195 0.000000000
mosaik196 0.000000000
mosaik197 0.000000000
mosaik198 0.000000000
mosaik199 0.000000000
mosaik200 0.000000000
mosaik201 0.000000000
mosaik202 0.000000000
mosaik203 0.000000000
mosaik204 0.000000000
mosaik205 0.000000000
mosaik206 0.000000000
mosaik207 0.000000000
mosaik208 0.000000000
mosaik209 0.000000000
mosaik210 0.000000000
mosaik211 0.000000000
mosaik212 0.000000000
mosaik213 0.000000000
mosaik214 0.000000000
mosaik215 0.000000000
mosaik216 0.000000000
mosaik217 0.000000000
mosaik218 0.000000000
mosaik219 0.000000000
mosaik220 0.000000000
mosaik221 0.000000000
mosaik222 0.000000000
mosaik223 0.000000000
mosaik224 0.000000000
mosaik225 0.000000000
mosaik226 0.000000000
mosaik227 0.000000000
mosaik228 0.000000000
mosaik229 0.000000000
mosaik230 0.000000000
mosaik231 0.000000000
mosaik232 0.000000000
mosaik233 0.000000000
mosaik234 0.015558794
mosaik235 0.000000000
mosaik236 0.000000000
mosaik237 0.000000000
mosaik238 0.000000000
mosaik239 0.000000000
mosaik240 0.000000000
mosaik241 0.000000000
mosaik242 0.000000000
mosaik243 0.000000000
mosaik244 0.000000000
mosaik245 0.000000000
mosaik246 0.000000000
mosaik247 0.000000000
mosaik248 0.000000000
mosaik249 0.000000000
mosaik250 0.000000000
mosaik251 0.000000000
mosaik252 0.000000000
mosaik253 0.000000000
mosaik254 0.000000000
mosaik255 0.000000000
mosaik256 0.000000000
mosaik257 0.000000000
mosaik258 0.000000000
mosaik259 0.000000000
mosaik260 0.000000000
mosaik261 0.000000000
mosaik262 0.000000000
mosaik263 0.000000000
mosaik264 0.000000000
mosaik265 0.000000000
mosaik266 0.000000000
mosaik267 0.000000000
mosaik268 0.000000000
mosaik269 0.000000000
mosaik270 0.000000000
mosaik271 0.000000000
mosaik272 0.000000000
mosaik273 0.000000000
mosaik274 0.000000000
mosaik275 0.000000000
mosaik276 0.000000000
mosaik277 0.044130256
mosaik278 0.000000000
mosaik279 0.000000000
mosaik280 -8.148026851
mosaik281 0.000000000
mosaik282 0.000000000
mosaik283 0.000000000
mosaik284 0.000000000
mosaik285 0.000000000
mosaik286 0.000000000
mosaik287 0.000000000
mosaik288 0.000000000
mosaik289 0.000000000
mosaik290 0.000000000
mosaik291 0.000000000
mosaik292 0.000000000
mosaik293 0.000000000
mosaik294 0.000000000
mosaik295 0.000000000
mosaik296 0.000000000
mosaik297 0.000000000
mosaik298 0.000000000
mosaik299 0.000000000
mosaik300 0.000000000
mosaik301 0.000000000
mosaik302 0.000000000
mosaik303 0.000000000
mosaik304 0.000000000
mosaik305 0.000000000
mosaik306 0.000000000
mosaik307 0.000000000
mosaik308 0.000000000
mosaik309 0.000000000
mosaik310 0.000000000
mosaik311 0.000000000
mosaik312 0.000000000
mosaik313 0.000000000
mosaik314 0.000000000
mosaik315 0.000000000
mosaik316 0.000000000
mosaik317 0.000000000
mosaik318 0.000000000
mosaik319 0.000000000
mosaik320 0.000000000
mosaik321 0.000000000
mosaik322 0.000000000
mosaik323 0.000000000
mosaik324 0.000000000
mosaik325 0.000000000
mosaik326 0.000000000
mosaik327 0.000000000
mosaik328 0.000000000
mosaik329 0.000000000
mosaik330 0.000000000
mosaik331 0.000000000
mosaik332 0.000000000
mosaik333 0.000000000
mosaik334 0.000000000
mosaik335 0.000000000
mosaik336 0.000000000
mosaik337 0.000000000
mosaik338 0.000000000
mosaik339 0.000000000
mosaik340 0.000000000
mosaik341 0.000000000
mosaik342 0.000000000
mosaik343 0.000000000
mosaik344 0.000000000
mosaik345 0.000000000
mosaik346 0.000000000
mosaik347 0.000000000
mosaik348 0.000000000
mosaik349 0.000000000
mosaik350 0.000000000
mosaik351 0.000000000
mosaik352 0.000000000
mosaik353 0.000000000
mosaik354 0.000000000
mosaik355 0.000000000
mosaik356 0.000000000
mosaik357 0.000000000
mosaik358 0.000000000
mosaik359 0.000000000
mosaik360 0.000000000
mosaik361 0.000000000
mosaik362 0.000000000
mosaik363 0.000000000
mosaik364 0.000000000
mosaik365 0.000000000
mosaik366 0.000000000
mosaik367 0.000000000
mosaik368 0.000000000
mosaik369 0.000000000
mosaik370 0.000000000
mosaik371 0.000000000
mosaik372 0.000000000
mosaik373 0.000000000
mosaik374 0.000000000
mosaik375 0.000000000
mosaik376 0.000000000
mosaik377 0.000000000
mosaik378 0.000000000
mosaik379 0.000000000
mosaik380 0.000000000
mosaik381 0.000000000
mosaik382 0.000000000
mosaik383 0.000000000
mosaik384 0.000000000
mosaik385 0.000000000
mosaik386 0.000000000
mosaik387 0.000000000
mosaik388 0.000000000
mosaik389 0.000000000
mosaik390 0.000000000
mosaik391 0.000000000
mosaik392 0.000000000
mosaik393 0.000000000
mosaik394 0.000000000
mosaik395 0.000000000
mosaik396 -35.261808063
mosaik397 0.000000000
mosaik398 0.000000000
mosaik399 0.000000000
mosaik400 0.000000000
mosaik401 0.000000000
mosaik402 0.000000000
mosaik403 0.000000000
mosaik404 0.000000000
mosaik405 0.000000000
mosaik406 0.000000000
mosaik407 0.000000000
mosaik408 0.000000000
mosaik409 0.000000000
mosaik410 0.000000000
mosaik411 0.000000000
mosaik412 0.000000000
mosaik413 0.000000000
mosaik414 0.000000000
mosaik415 0.000000000
mosaik416 0.000000000
mosaik417 0.000000000
mosaik418 0.000000000
mosaik419 0.000000000
mosaik420 0.000000000
mosaik421 0.000000000
mosaik422 0.000000000
mosaik423 0.000000000
mosaik424 0.000000000
mosaik425 0.000000000
mosaik426 0.000000000
mosaik427 0.000000000
mosaik428 0.000000000
mosaik429 0.000000000
mosaik430 0.000000000
mosaik431 0.000000000
mosaik432 0.000000000
mosaik433 0.000000000
mosaik434 0.000000000
mosaik435 0.000000000
mosaik436 0.000000000
mosaik437 0.000000000
mosaik438 0.000000000
mosaik439 0.000000000
mosaik440 0.000000000
mosaik441 0.000000000
mosaik442 0.000000000
mosaik443 0.000000000
mosaik444 0.000000000
mosaik445 0.000000000
mosaik446 0.000000000
mosaik447 0.000000000
mosaik448 0.000000000
mosaik449 0.000000000
mosaik450 0.000000000
mosaik451 0.000000000
mosaik452 0.000000000
mosaik453 0.000000000
mosaik454 0.000000000
mosaik455 0.000000000
mosaik456 0.000000000
mosaik457 0.000000000
mosaik458 0.000000000
mosaik459 55.635488670
mosaik460 0.000000000
mosaik461 0.000000000
mosaik462 0.000000000
mosaik463 0.000000000
mosaik464 0.000000000
mosaik465 0.000000000
mosaik466 0.000000000
mosaik467 0.000000000
mosaik468 0.000000000
mosaik469 0.000000000
mosaik470 0.000000000
mosaik471 0.000000000
mosaik472 0.000000000
mosaik473 0.000000000
mosaik474 0.000000000
mosaik475 0.000000000
mosaik476 0.000000000
mosaik477 0.000000000
mosaik478 0.000000000
mosaik479 0.000000000
mosaik480 0.000000000
mosaik481 0.000000000
mosaik482 0.000000000
mosaik483 0.000000000
mosaik484 0.000000000
mosaik485 0.000000000
mosaik486 0.000000000
mosaik487 0.000000000
mosaik488 0.000000000
mosaik489 0.000000000
mosaik490 0.000000000
mosaik491 0.000000000
mosaik492 0.000000000
mosaik493 0.000000000
mosaik494 0.000000000
mosaik495 0.000000000
mosaik496 0.000000000
mosaik497 0.000000000
mosaik498 0.000000000
mosaik499 0.000000000
mosaik500 0.000000000What we want: the non-zero variable names!
- Getting the names of the variables is more complicated than it should be
Code
s1 variable
(Intercept) 0.5686581 (Intercept)
mosaik1 0.0000000 mosaik1
mosaik2 0.0000000 mosaik2
mosaik3 0.0000000 mosaik3
mosaik4 0.0000000 mosaik4
mosaik5 0.0000000 mosaik5Code
[1] "(Intercept)" "mosaik39" "mosaik234" "mosaik277" "mosaik280"
[6] "mosaik396" "mosaik459" One more step: remove the Intercept!
- We don’t want the name of the intercept
- All of the packages we use will add that automatically
How do we use this with ebp?
- In EBP, we need a
formula - How do we turn this into a formula?
- We need to add the outcome variable (poor)
andcombine the predictors with+
- We need to add the outcome variable (poor)
Finally: estimating the model
Code
library(povmap) # I like to use povmap instead of emdi (personal preference)
# get "area" variable
predictors$TA_CODE <- substr(predictors$EA_CODE, 1, 5)
data$TA_CODE <- substr(data$EA_CODE, 1, 5)
ebp <- ebp(fixed = ebpformula, # the formula
pop_data = predictors, # the population data
pop_domains = "TA_CODE", # the domain (area) name in the population data
smp_data = data, # the sample data
smp_domains = "TA_CODE", # the domain (area) name in the sample data
transformation = "arcsin", # I'm going to use the arcsin transformation
weights = "total_weights", # sample weights
weights_type = "nlme", # weights type
MSE = TRUE, # variance? yes please
L = 0) # this is a new thing in povmap: "analytical" point estimates. much faster!Time difference of 1.12 secsEmpirical Best Prediction
Call:
ebp(fixed = poor ~ mosaik39 + mosaik234 + mosaik277 + mosaik280 +
mosaik396 + mosaik459, pop_data = predictors, pop_domains = "TA_CODE",
smp_data = data, smp_domains = "TA_CODE", L = 0, transformation = "arcsin",
MSE = TRUE, weights = "total_weights", weights_type = "nlme")
Out-of-sample domains: 27
In-sample domains: 49
Out-of-sample subdomains: 0
In-sample subdomains: 0
Sample sizes:
Units in sample: 107
Units in population: 2911
Min. 1st Qu. Median Mean 3rd Qu. Max.
Sample_domains 1 1.0 2 2.183673 3 7
Population_domains 1 5.5 18 38.302632 44 300
Explanatory measures for the mixed model:
Marginal_R2 Conditional_R2 Marginal_Area_R2 Conditional_Area_R2
0.2552368 0.4115571 0.4196524 0.6605383
Residual diagnostics for the mixed model:
Skewness Kurtosis Shapiro_W Shapiro_p
Error -0.5081831 3.993441 0.9761079 0.05063839
Random_effect 0.1489666 4.322010 0.9666000 0.17669350
Estimated variance of random effects:
Variance
Error 0.028541022
Random_effect 0.007581942
ICC: 0.2098926
Shrinkage factors
gamma.Min. gamma.1st.Qu. gamma.Median gamma.Mean gamma.3rd.Qu.
TA_CODE 0.2098926 0.2098926 0.3272247 0.3240028 0.4083161
gamma.Max.
TA_CODE 0.6419322
Transformation:
Transformation Shift_parameter
arcsin 0Some results
Some results
Empirical Best Prediction
Call:
ebp(fixed = poor ~ mosaik39 + mosaik234 + mosaik277 + mosaik280 +
mosaik396 + mosaik459, pop_data = predictors, pop_domains = "TA_CODE",
smp_data = data, smp_domains = "TA_CODE", L = 0, transformation = "arcsin",
MSE = TRUE, weights = "total_weights", weights_type = "nlme")
Out-of-sample domains: 27
In-sample domains: 49
Out-of-sample subdomains: 0
In-sample subdomains: 0
Sample sizes:
Units in sample: 107
Units in population: 2911
Min. 1st Qu. Median Mean 3rd Qu. Max.
Sample_domains 1 1.0 2 2.183673 3 7
Population_domains 1 5.5 18 38.302632 44 300
Explanatory measures for the mixed model:
Marginal_R2 Conditional_R2 Marginal_Area_R2 Conditional_Area_R2
0.2552368 0.4115571 0.4196524 0.6605383
Residual diagnostics for the mixed model:
Skewness Kurtosis Shapiro_W Shapiro_p
Error -0.5081831 3.993441 0.9761079 0.05063839
Random_effect 0.1489666 4.322010 0.9666000 0.17669350
Estimated variance of random effects:
Variance
Error 0.028541022
Random_effect 0.007581942
ICC: 0.2098926
Shrinkage factors
gamma.Min. gamma.1st.Qu. gamma.Median gamma.Mean gamma.3rd.Qu.
TA_CODE 0.2098926 0.2098926 0.3272247 0.3240028 0.4083161
gamma.Max.
TA_CODE 0.6419322
Transformation:
Transformation Shift_parameter
arcsin 0Some results
Indicator/s: Mean
Domain Mean Mean_MSE Mean_CV
1 10101 0.3669687 0.002831292 0.1449984
2 10102 0.3526612 0.004281904 0.1855499
3 10103 0.2829257 0.005202843 0.2549458
4 10104 0.3255265 0.010839984 0.3198365
5 10105 0.3564514 0.005060076 0.1995621
6 10106 0.3747787 0.004249448 0.1739367
7 10107 0.3518191 0.006747908 0.2334883
8 10108 0.3545512 0.011015621 0.2960230
9 10109 0.3971006 0.005494159 0.1866595
10 10110 0.2786403 0.029178299 0.6130360
11 10120 0.3480978 0.010632100 0.2962158
12 10201 0.3514189 0.003631778 0.1714883
13 10202 0.3344865 0.005765226 0.2270019
14 10203 0.4354745 0.003538805 0.1366047
15 10204 0.3489391 0.003373777 0.1664595
16 10205 0.3882951 0.004288408 0.1686499
17 10206 0.2853691 0.028118716 0.5876122
18 10220 0.1969776 0.003507940 0.3006832
19 10301 0.2713425 0.004896562 0.2578860
20 10302 0.2897256 0.004537121 0.2324895
21 10303 0.3130857 0.006657093 0.2606028
22 10304 0.2778572 0.004797752 0.2492857
23 10305 0.2522589 0.005147865 0.2844246
24 10306 0.2242030 0.003030354 0.2455304
25 10307 0.2762336 0.007211858 0.3074305
26 10308 0.2357340 0.005214202 0.3063174
27 10309 0.3058134 0.004771198 0.2258693
28 10310 0.2791358 0.007478376 0.3098047
29 10311 0.3226810 0.023091996 0.4709310
30 10312 0.1930645 0.004107526 0.3319614
31 10314 0.3159857 0.018484925 0.4302703
32 10320 0.2135503 0.006209553 0.3690030
33 10401 0.2920620 0.002668281 0.1768645
34 10402 0.2809861 0.012008835 0.3900009
35 10403 0.2862341 0.004815276 0.2424316
36 10404 0.3170657 0.004813032 0.2188065
37 10405 0.4033321 0.005354227 0.1814202
38 10406 0.2999634 0.017052315 0.4353348
39 10407 0.3129457 0.006060281 0.2487582
40 10408 0.3671935 0.006664314 0.2223221
41 10409 0.3303997 0.008374931 0.2769816
42 10410 0.2798613 0.006156521 0.2803656
43 10411 0.3233065 0.011491456 0.3315684
44 10412 0.3241912 0.016624524 0.3977164
45 10413 0.2911344 0.039942797 0.6864765
46 10420 0.1852033 0.003687328 0.3278742
47 10501 0.4279057 0.005303207 0.1701850
48 10502 0.3145523 0.004670334 0.2172607
49 10503 0.3572890 0.006584669 0.2271158
50 10504 0.2799923 0.009966906 0.3565613
51 10505 0.4006467 0.005227938 0.1804695
52 10506 0.3835653 0.005840463 0.1992436
53 10507 0.3457535 0.011226908 0.3064527
54 10508 0.4131363 0.005102155 0.1728954
55 10509 0.4611727 0.012754984 0.2448930
56 10510 0.4246205 0.007247554 0.2004909
57 10511 0.3551883 0.007466637 0.2432785
58 10512 0.3095397 0.032565027 0.5829876
59 10520 0.2261618 0.013804918 0.5195146
60 10601 0.2232001 0.009288918 0.4318055
61 10620 0.2148553 0.038237449 0.9101196
62 10701 0.1332353 0.019293466 1.0425232
63 10702 0.1960373 0.008067828 0.4581838
64 10703 0.1714235 0.020416275 0.8335235
65 10704 0.1271540 0.020051093 1.1136247
66 10705 0.2342148 0.012515007 0.4776405
67 10706 0.1794717 0.010035269 0.5581726
68 10707 0.1332409 0.005456031 0.5543713
69 10708 0.1786822 0.006756107 0.4600098
70 10709 0.1986310 0.014408767 0.6043193
71 10710 0.2180149 0.004968069 0.3233015
72 10711 0.1896444 0.009930312 0.5254622
73 10712 0.1853390 0.029157502 0.9213156
74 10713 0.2527483 0.009113402 0.3777044
75 10714 0.1688222 0.010341302 0.6023628
76 10715 0.1517842 0.027500987 1.0925659Let’s put those results on a map!
Let’s put those results on a map!
What does this look like if we DON’T transform?
Code
ebpnotransform <- ebp(fixed = ebpformula, # the formula
pop_data = predictors, # the population data
pop_domains = "TA_CODE", # the domain (area) name in the population data
smp_data = data, # the sample data
smp_domains = "TA_CODE", # the domain (area) name in the sample data
transformation = "no",
weights = "total_weights", # sample weights
weights_type = "nlme", # weights type
MSE = TRUE, # variance? yes please
L = 0) # this is a new thing in povmap: "analytical" point estimates. much faster!Time difference of 1.11 secsThe residuals
Comparing the models
No transformation
Empirical Best Prediction
Call:
ebp(fixed = poor ~ mosaik39 + mosaik234 + mosaik277 + mosaik280 +
mosaik396 + mosaik459, pop_data = predictors, pop_domains = "TA_CODE",
smp_data = data, smp_domains = "TA_CODE", L = 0, transformation = "no",
MSE = TRUE, weights = "total_weights", weights_type = "nlme")
Out-of-sample domains: 27
In-sample domains: 49
Out-of-sample subdomains: 0
In-sample subdomains: 0
Sample sizes:
Units in sample: 107
Units in population: 2911
Min. 1st Qu. Median Mean 3rd Qu. Max.
Sample_domains 1 1.0 2 2.183673 3 7
Population_domains 1 5.5 18 38.302632 44 300
Explanatory measures for the mixed model:
Marginal_R2 Conditional_R2 Marginal_Area_R2 Conditional_Area_R2
0.2140398 0.4022848 0.363542 0.6916272
Residual diagnostics for the mixed model:
Skewness Kurtosis Shapiro_W Shapiro_p
Error 0.08332729 3.198589 0.9944134 0.94485679
Random_effect 0.61867234 4.133345 0.9540069 0.05375011
Estimated variance of random effects:
Variance
Error 0.018896174
Random_effect 0.005951178
ICC: 0.2395095
Shrinkage factors
gamma.Min. gamma.1st.Qu. gamma.Median gamma.Mean gamma.3rd.Qu.
TA_CODE 0.2395095 0.2395095 0.3657341 0.35956 0.4499858
gamma.Max.
TA_CODE 0.6800417
Transformation: No transformation Arcsin transformation
Empirical Best Prediction
Call:
ebp(fixed = poor ~ mosaik39 + mosaik234 + mosaik277 + mosaik280 +
mosaik396 + mosaik459, pop_data = predictors, pop_domains = "TA_CODE",
smp_data = data, smp_domains = "TA_CODE", L = 0, transformation = "arcsin",
MSE = TRUE, weights = "total_weights", weights_type = "nlme")
Out-of-sample domains: 27
In-sample domains: 49
Out-of-sample subdomains: 0
In-sample subdomains: 0
Sample sizes:
Units in sample: 107
Units in population: 2911
Min. 1st Qu. Median Mean 3rd Qu. Max.
Sample_domains 1 1.0 2 2.183673 3 7
Population_domains 1 5.5 18 38.302632 44 300
Explanatory measures for the mixed model:
Marginal_R2 Conditional_R2 Marginal_Area_R2 Conditional_Area_R2
0.2552368 0.4115571 0.4196524 0.6605383
Residual diagnostics for the mixed model:
Skewness Kurtosis Shapiro_W Shapiro_p
Error -0.5081831 3.993441 0.9761079 0.05063839
Random_effect 0.1489666 4.322010 0.9666000 0.17669350
Estimated variance of random effects:
Variance
Error 0.028541022
Random_effect 0.007581942
ICC: 0.2098926
Shrinkage factors
gamma.Min. gamma.1st.Qu. gamma.Median gamma.Mean gamma.3rd.Qu.
TA_CODE 0.2098926 0.2098926 0.3272247 0.3240028 0.4083161
gamma.Max.
TA_CODE 0.6419322
Transformation:
Transformation Shift_parameter
arcsin 0One more thing: we can go higher!
Code
predictors$dist <- substr(predictors$TA_CODE, 1, 3)
data$dist <- substr(data$TA_CODE, 1, 3)
ebpnotransform <- ebp(fixed = ebpformula, # the formula
pop_data = predictors, # the population data
pop_domains = "TA_CODE", # the domain (area) name in the population data
aggregate_to = "dist", # aggregate to district!
smp_data = data, # the sample data
smp_domains = "TA_CODE", # the domain (area) name in the sample data
transformation = "no",
weights = "total_weights", # sample weights
weights_type = "nlme", # weights type
MSE = TRUE, # variance? yes please
L = 0) # this is a new thing in povmap: "analytical" point estimates. much faster!Time difference of 1.02 secs